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Abstract;  This  paper  introduces  a  time-dependent  model  of  regional  barrier 
island  breaching  for  multiple  openings,  whether  as  permanent  inlets  or  as 
new  breaches  to  the  same  bay.  The  model  allows  an  arbitrary  number  of 
breaches  and  is  forced  by  tide,  storm  surge,  and  wave  set  up.  Limitations  on 
the  hydrodynamics,  hence  calculated  breach  evolution,  are  those  associated 
with  the  Keulegan  inlet  model.  The  model  includes  possible  closure  or 
limited  breaching  by  shoaling  of  the  breach  channel  by  longshore  transport. 
Capabilities  of  the  multiple  breaching  model  are  examined  through  eight 
sensitivity  tests. 

INTRODUCTION 

A  breach  in  a  barrier  island  occurs  if  inundation  with  strong  flow  develops  a  channel 
through  the  barrier.  The  inundation  may  be  forced  from  the  ocean  side  or  the  bay  side. 
During  a  breach,  the  initial  channel  typically  expands  vertically  and  horizontally  in 
exponential  manner  (Kraus  2003)  to  dynamic  equilibrium  dimensions.  Forcing  of  the 
water  flow  in  the  breach  is  caused  by  a  difference  between  bay  and  ocean  water  levels, 
influenced  by  wave  runup  and  overwash. 

High  ocean  water  surface  elevations  and  large  waves  associated  with  severe  storms  can 
produce  multiple  breaches  along  a  coastal  barrier  island  complex  at  regional  spatial 
scales.  Occurrences  of  multiple  breach  channels  cut  by  storms  has  been  reported  (e.g., 
Price  1947;  Terchunian  and  Merkert  1995;  Wamsley  and  Hathaway  2004;  Wutkowski 
2004;  Freeman  et  al.  2004),  but  there  remains  a  paucity  of  quantitative  barrier  island 
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breach  field  data.  Reviews  of  coastal  barrier  island  breaching  causes,  processes,  and 
remediation  are  given  by  Kraus  and  Wamsley  (2003)  and  Wamsley  and  Kraus  (2005). 

Kraus  and  Hayashi  (2005)  developed  a  process-based  numerical  model  to  simulate 
morphologic  evolution  of  coastal  barrier  island  breaching  for  engineering  use.  The 
objective  of  the  present  study  is  to  extend  the  capabilities  of  the  Kraus  and  Hayashi 
model  to  include  development  and  evolution  of  multiple  breaches  in  areas  along  a  barrier 
island,  which  may  include  presence  of  permanent  inlets.  This  multiple-breach 
morphological  model  is  applied  on  regional  spatial  scales  by  integration  with  the 
Cascade  model  of  regional  sediment  transport  and  shoreline  change  (Larson  et  al.  2002; 
Larson  and  Kraus  2003;  Larson  et  al.  2006;  Connell  and  Kraus  2006;  Larson  et  al.  2007). 

COASTAL  BREACH  MODEL 

The  morphologic  breaching  model  is  based  on  analytical  models  for  inlet  stability  and 
inlet  channel  cross-sectional  area  (Kraus  1998,  2003),  and  it  has  been  extended  to 
include  morphologic  evolution  of  a  rectangular  barrier  island  (Kraus  2003)  and  a  stacked 
rectangle  pyramidal  barrier  island  (Kraus  and  Hayashi  2005).  Studies  have  examined 
predictions  of  the  coastal  breach  morphology  model  by  comparison  to  field  and 
laboratory  breaches  (Kraus  and  Hayashi  2005;  Wamsley  et  al.  2006;  Wamsley  et  al. 
2007).  The  coastal  breaching  model  simulates  breaching  by  iteratively  solving  the  1-D 
Keulegan  (1967)  equations  for  inlet  hydrodynamics  to  determine  breach  current  velocity 
and  bay  elevation  through  time  based  on  breach  channel  cross-sectional  area  and  area  of 
the  bay.  The  coastal  breaching  model  calculates  under  the  assumptions  that  the  bay  must 
accommodate  uniform  water  surface  deviation,  the  bay  is  composed  of  vertical  walls  so 
that  volume  and  surface  area  are  directly  correlated,  and  the  breaches  do  not  overlap. 

If  multiple  breach  channels  are  calculated  to  occur  in  the  simulation,  the  Keulegan 
equations  become: 
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for  momentum,  and: 

P) 

z=i  cit 

for  conservation  of  water  volume,  where  Ui  =  depth-averaged  and  inlet- length  integrated 
current  velocity  at  breach  or  inlet  i  of  N  openings;  t  =  time;  g  =  gravitational 
acceleration;  L  =  width  of  breach  through  the  barrier  island;  r|Q  =  water  surface 

deviations  from  Mean  Sea  Level  (MSL)  in  the  ocean;  Ken  and  Kex  =  entrance  and  exit 
losses,  respectively  (non-dimensional);  Rh= hydraulic  radius  of  the  breach  or  inlet  (m); 
Ac = breach  or  inlet  channel  cross-sectional  area  below  MSL  (m  );  ABay  =  surface  area  of 
the  bay  (m2).  The  bottom  friction  coefficient  c/is  represented  by  Mannings  formula, 
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in  which  n  =  Mannings  coefficient  typically  set  as  0.025  s/m  ;  and  h  =  water  depth. 

Sediment  transport  through  the  breach  is  calculated  employing  the  Watanabe  et  al. 
(1991)  total  load  formula: 


qB  -  a—  T—U  (4) 

Pg 

in  which  a  =  empirical  coefficient  typically  of  order  0.1  to  10  (non-dimensional); 
Tm  =  cfp[\u\u),  the  time-averaged  bottom  shear  stress  for  waves  and  current; 

rc  =  (ps  -  p)gd50x ¥c ,  the  critical  shear  for  sediment  motion,  where  ps  is  the  density  of 
sand  (kg/m3),  p  is  the  density  of  water  (kg/m3),  d 50  =  median  grain  size  diameter  (m); 
and  'hc  =  critical  Shields  parameter,  set  to  0.05  for  fine  to  medium  sand. 


If  waves  are  not  included  in  the  breaching  simulation,  Eq.  4  becomes: 
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Sediment  transport  along  the  walls  of  the  channel  is  calculated  as  a  fraction  of  the 
transport  along  the  breach  floor.  A  typical  value  for  side  transport  is  approximately  70% 
of  the  transport  along  the  channel. 


qs  =  fjch  (6) 

where  p  is  the  calibration  factor  for  the  channel  wall  fraction  (e.g.,  0.7). 

The  wave  setup  contribution  to  mean  water  level  at  the  shoreline  is  calculated  following 
Stockdon  et  al.  (2006): 


(7) 

where  Pf  =  beach  foreshore  slope;  H0  =  deep-water  wave  height  (m);  L0  =  deep-water 

wavelength  (m).  Kraus  and  Hayashi  2005;  Wamsley  et  al.  (2006)  discuss  calculation  of 
channel  infdling  and  shoaling  due  to  longshore  transport  estimates  calculated  from  wave 
input  data.  Longshore  transport  is  implemented  in  the  coastal  breaching  model  with  the 
CERC  formula  in  the  Coastal  Engineering  Manual  (Headquarters,  U.S.  Army  Corps  of 
Engineers  2002): 
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where  K  =  dimensionless  sediment  transport  coefficient  based  on  wave  energy, 
K=breaker  index  (non-dimensional  ~  0.8),  n  =  sediment  porosity,  Hb  rms  =  root-mean- 

square  wave  height  at  breaking  point  (m),  and  ab  =  wave  breaking  angle. 


The  present  study  extends  the  capabilities  of  the  morphologic  breaching  model  by 
implementing  multiple  breaches  with  variable  (in  time  and  space)  breach  forcing. 
Furthermore,  breaching  model  connectivity  with  the  Cascade  model  is  now  possible  by 
employing  a  subset  of  the  wave  and  water  surface  elevation  data  used  as  inputs  to 
Cascade,  and  spatially  and  temporally  variable  longshore  sediment  transport.  Cascade 
output  can  control  breach  infilling  and  shoaling  rates  by  longshore  transport  in  place  of 
Eq.  8  in  the  breaching  model,  subject  of  a  planned  future  paper. 

METHOD 

An  idealized  alongshore  barrier  dune-crest  topography  (Fig.  1)  was  developed  to 
simulate  breaching  at  topographic  lows  in  the  dunes.  This  topography  can  be  arbitrary, 
and  the  topographic  lows  function  as  potential  pilot  channels  for  the  breaching  model. 
Pilot  channels  of  different  elevations  may  exist,  providing  selective  inundation  for 
different  water  levels  at  various  times  in  a  storm  and  thereafter.  Morphologic  evolution 
of  breach  channels  occurs  only  if  the  pilot  channel  is  inundated.  Local  forcing,  including 
longshore  transport,  waves,  water  level,  and  shoreline  position  for  the  breaches  is 
provided  at  regional  scale  by  the  Cascade  model.  Sensitivity  tests  are  conducted  to 
isolate  some  of  the  complex  processes  that  force  change  in  breach  morphology  and  to 
evaluate  the  functionality  of  the  model  under  hypothetical  forcing. 
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Eight  14-day  simulations  of  the  evolution  of  multiple  breach  channels  were  conducted 
(Table  1)  with  the  idealized  barrier  alongshore  topography.  Dune  crest  elevation  was 
3.5  m,  and  four  different  pilot  channels  10  m  wide  by  50  m  long  were  cut  with  elevations 
of  3,  2.5,  2,  and  1.5  m.  The  location  of  the  pilot  channels  along  the  barrier  island  is 
arbitrary  and  influenced  only  where  input  forcing  is  (optionally)  extracted  from  Cascade. 
It  is  assumed  that  the  breaches  cannot  converge,  and  therefore,  pilot  channels  cannot  be 
located  at  points  along  the  barrier  where  the  potential  exists  for  two  or  more  breaches  to 
merge  into  one.  Bay  surface  area  was  held  constant  at  5xl08  m2.  The  time-step  for  the 
simulations  was  set  at  10  s.  Empirical  coefficients  for  sediment  transport  along  the 
bottom  and  along  the  side  of  the  breach  were  set  to  0.5  and  0.7,  respectively.  Sediment 
porosity  was  0.4,  and  the  //-coefficient  of  the  CERC  formula  was  set  to  0.39  (based  on 
significant  wave  height).  Finally,  water  (sediment)  density  was  set  to  1,025  kg/m3 
(2,650  kg/m3). 

The  first  sensitivity  test  (Case  HD)  held  a  constant  0.2-m  head  difference  along  the  entire 
length  of  the  barrier  generated  by  a  constant  ocean  water  surface  elevation  at  3-m  above 
MSL  and  a  constant  bay  water  surface  elevation  at  2.8-m  above  MSL.  Two  other  cases 
were  conducted  employing  the  same  constant  head  difference  water  level  forcing 
between  ocean  and  bay,  added  wave  forcing  to  simulate  wave  setup  forcing  and  sediment 
infilling  due  to  longshore  transport.  Case  HD-H1T8  D0  includes  waves  with  constant 
Ho=  1  m,  To  =  8  s,  and  a0  =  0  deg,  and  Case  HD-H1T8D30  includes  waves  with 

constant  Ho  =  1  in,  Tq=  8  s,  and  a0  =  30  deg. 

Three  sensitivity  tests  were  conducted  forced  with  a  1-m  amplitude  M2  (12. 42 -hr  period) 
tidal  wave  with  90-deg  phase  lag  between  the  ocean  and  bay  levels,  applied  along  the 
entire  length  of  the  barrier  with  a  vertical  datum  shift  of  2-m  above  MSL.  These  cases 
are  listed  in  Table  1  as  Case  M2,  Case  M2-H1T8D0,  and  Case  M2-H1T8D30.  Two 
cases  also  included  wave  forcing  to  estimate  wave  setup  and  sediment  infilling  due  to 
longshore  transport.  Case  M2-H1T8  D0  includes  waves  with  constant  H{)  =  1  m,  7o  =  8 
s,  and  a0=  0  deg,  and  Case  M2-H1_T8_D30  includes  waves  with  constant  Hq=  1  m, 

To  =  8  s,  and  a0  =  30  deg. 

The  two  remaining  sensitivity  tests  examined  response  of  the  model  to  a  synthetic  storm 
surge.  The  first  of  the  storm  surge  simulations  (Case  Surge  in  Table  1)  was  forced  with 
a  constant  water  level  1  m  above  MSL  until  day  5,  after  which  a  2  m  simulated  parabolic 
storm  surge  extending  over  2  days  was  added  to  initiate  breach  inundation  (Fig.  2).  This 
surge  was  set  out  of  phase  for  each  breach  channel  so  that  each  breach  would  be  initiated 
6  hr  after  the  previous  breach.  The  second  surge  case  (Case  Surge+M2  in  Table  1) 
consisted  of  the  same  synthetic  surge  as  Case  Surge,  but  the  surge  arrival  was  equal  for 
all  breaches,  and  a  simple  1  -m  M2  tide  with  90-deg  phase  lag  between  the  ocean  and  bay 
levels  was  added  to  the  surge  (Fig.  3).  This  M2  tide  differed  from  Case  M2  in  that  it  did 
not  reach  2  m  above  MSL  to  prevent  breach  inundation  prior  to  the  surge  arrival. 
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Fig.  2.  Ocean  and  bay  water  level  -  Case 
Surge 


Fig.  3.  Ocean  and  bay  water  level  -  Case 
Surge+M2 


Table  1.  Regional  Breach  Model  Sensitivity  Test  Cases 


Case  Name 

Water  level  forcing 

Deep  Water 
Wave  Height 
(Ho),  m 

Deep  Water 
Wave  Period 
(To),  s 

Deep  Water 
Wave  Deg 
(«o),  deg 

HD 

Constant  ocean  =  3  m 
constant  bay  =  2.8  m 

N/A 

N/A 

N/A 

HD-H1_T8_D0 

Constant  ocean  =  3  m 
constant  bay  =  2.8  m 

1 

8 

0 

HD-H1_T8_D30 

Constant  ocean  =  3  m 
constant  bay  =  2.8  m 

1 

8 

30 

M2 

M2-tide  amp:  1  m 
phase:  ocean  =  0° 
phase:  bay  =  90° 

N/A 

N/A 

N/A 

M2-H1_T8_D0 

M2-tide  amp:  1  m 
phase:  ocean  =  0° 
phase:  bay  =  90° 

1 

8 

0 

M2-H1_T8_D30 

M2-tide  amp:  1  m 
phase:  ocean  =  0° 
phase:  bay  =  90° 

1 

8 

30 

Surge 

2-m  storm  surge 
w/6-min  arrival  delay 
for  consecutive  breach 

N/A 

N/A 

N/A 

Surge+M2 

2-m  storm  surge  + 
M2-tide  amp:  1-m 
phase:ocean  =  0° 
phase:  bay  =  90° 

N/A 

N/A 

N/A 
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RESULTS 

Calculated  evolution  of  breach  channel  depth,  width,  total  transport  through  the  four 
breaches,  and  cross-sectional  area  are  shown  in  Figs.  4-7.  Breaches  that  started  in  pilot 
channels  at  lower  elevation  (e.g.,  elevation  z  =  1.5  m,  z  =  2  m)  were  cut  deeper  (Fig.  4) 
and  at  a  faster  rate  than  those  that  developed  in  pilot  channels  at  higher  elevations  (e.g., 
z  =  2.5  m,  z  =  3  m).  Conversely,  the  lateral  growth  of  higher  elevation  channels  in 
Fig.  5  is  greater  than  the  deeper,  low-elevation  channels.  The  greater  lateral  growth  is 
explained  by  increased  shear  stresses  by  bottom  friction  in  the  sediment  transport  and 
hydrodynamic  models.  Breaching  is  not  induced  at  the  pilot  channel  stating  at  z  =  3  m, 
because  this  elevation  is  the  same  as  the  ocean  water  surface  elevation.  Stair-shaped 
behavior  observed  in  Fig.  6  is  an  artifact  of  the  layered  pyramidal  barrier  island, 
corresponding  to  the  transport  rates  being  reduced  as  a  result  of  channel  cuts  into 
incrementally  longer  (in  cross-shore  orientation)  barrier  island  layers.  As  breach  channel 
depth  increases,  the  transport  can  increase  in  magnitude  and  duration. 

Figure  7  shows  the  breach  cross-sectional  area  over  the  simulation  period  for  Case  HD. 
Breach  channel  cross-sectional  area  is  a  representation  of  the  channel  dimensions  and 
serves  as  a  good  indicator  of  channel  stability  and  channel  evolution  morphology.  A 
delayed  and  extended  incipient  breaching  stage  occurs  in  the  breach  channel  initially  at 
z  =  2.5  m  relative  to  the  channels  at  z  =  2  and  z  =  1 .5  m.  This  is,  again,  explained  by  the 
increased  shear  stress  associated  with  shallow  water  over  the  channel  bottom.  A 
shallower  channel  is  forced  to  grow  laterally  if  sediment  transport  rates  (Fig.  6)  increase 
due  to  greater  shear  stress  applied  over  a  larger  channel  bottom  surface  area  than 
generated  in  deeper  breach  channels. 

Figures  8  and  9  compare  breach  channel  cross-sectional  area  for  Cases  HD-H1T8D0 
and  HD-H1T8D30,  respectively.  Figure  8  retains  many  of  the  same  patterns  observed 
in  the  cross-sectional  area  plots  of  Case  HD  (Fig.  7).  However,  because  of  the  presence 
of  waves  and  the  associated  wave  setup,  a  shorter  incipient  breaching  stage  is  observed 
at  the  breach  of  z  =  2.5  m,  and  the  channel  cross-sectional  area  for  each  breach  in  Fig.  8 
at  the  end  of  the  simulation  is  greater  than  those  observed  in  Fig.  7.  Figure  9  shows  an 
entirely  different  pattern.  In  fact,  due  to  the  large  incident  wave  angle  relative  to 
shoreline,  significant  longshore  transport  causes  most  of  the  channels  to  shoal.  Breach 
closure  does  not  occur  for  the  breach  at  z  =  1.5  m  because  this  breach  pilot  channel  is 
sufficiently  deep  to  avoid  closure  during  the  incipient  phase  of  the  breach. 
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Fig.  4.  Calculated  breach  depth  ,  Case  HD 


Fig.  5.  Calculated  breach  width,  Case  HD 
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Fig.  6.  Calculated  sediment  transport  rate  Fig.  7.  Calculated  Breach  channel  cross- 
through  the  four  breaches,  Case  HD  sectional  area  evolution,  Case  HD 


Fig.  8.  Calculated  Breach  channel  cross-  Fig.  9.  Calculated  Breach  channel  cross- 
sectional  area  evolution,  HD_H1_T8_D0  sectional  area  evolution,  HD_H1_T8_D30 


Figures  10-13  illustrate  breach  transport  and  cross-sectional  area  for  the  M2, 
M2-H1T8  D0,  and  M2-H1_T8_D30  cases.  Figure  10  demonstrates  the  periodicity  of 
sediment  transport  associated  with  the  M2  water  elevation  forcing.  Depending  on  level 
of  the  breach  channel  relative  to  the  forcing  water  level,  transport  rates  may  never 
approach  zero,  whereas  sometimes  the  breach  channel  may  be  completely  depleted  of 
water  as  seen  in  Fig.  10  (z  =  2.5  m).  This  switching  on  and  off  of  transport  rates 
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significantly  reduces  total  transport  through  the  breach  as  shown  by  the  relative 
difference  in  cross-sectional  area  plots  (Fig.  11)  between  z  =  2.5  m  and  the  much  greater 
cross-sectional  area  represented  in  the  z  =  1.5  m  and  z  =  2  m  pilot  channels.  Case 
M2-H1T8  D0  (Fig.  12)  also  demonstrates  similar  morphologic  response  as  Case  M2 
(Fig.  1 1),  but  with  slightly  greater  cross-sectional  area  and  shorter  incipient  breach  phase 
times.  It  appears  that,  in  this  case,  the  addition  of  wave  setup  has  increased  flow  through 
the  breaches,  thereby  increasing  the  channel  area.  Significant  differences  in  cross- 
sectional  area  are  achieved  in  Case  M2-H1T8D30  with  waves  approaching  at  an 
incident  angle  (Fig.  13).  This  case  suggests  that  breach  channels  can  be  overwhelmed  by 
longshore  sediment  transport,  and  infilling  occurs  to  the  point  of  breach  closure  where 
breach  evolution  is  constant,  and  the  channels  are  no  longer  inundated. 
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Fig.  10.  Calculated  sediment  transport  rate  through  the  four  breaches,  Case  M2 


Fig.  11.  Calculated  breach  channel  cross-sectional  area  evolution,  Case  M2 
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Fig.  12.  Breach  channel  cross-sectional  area  evolution,  Case  M2-FI1_T8_D0 


Fig.  13.  Breach  channel  cross-sectional  area  evolution,  Case  M2-FI1_T8_D30 

Figures  14-17  illustrate  breach  transport  and  cross-sectional  area  for  the  Surge  and  the 
Surge+M2  cases.  Because  of  the  rapid  increase  in  water  level  from  the  simulated  storm 
surge,  peaks  in  sediment  transport  (Fig.  14)  are  present  immediately  after  inundation. 
These  peaks  grow  progressively  as  pilot  channels  of  initial  elevations  z  =  1 .5  m,  z  =  2  m, 
and  z  =  2.5  m  are  separately  inundated  and  expanded.  The  growth  of  the  peaks 
corresponds  to  the  increase  in  incipient  transport  according  to  the  depth  of  water  in  the 
breach  channel.  Shallow  channels  can  transport  more  sediment  because  of  larger  bottom 
shear  stress.  Figure  15  indicates  a  rapid  increase  in  cross-sectional  area  associated  with 
the  time  of  the  storm  surge. 
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Figure  16  illustrates  the  difference  caused  by  addition  of  a  tide  to  the  storm  surge. 
Because  of  the  higher  frequency  M2  tide  relative  to  the  surge,  sediment  transport 
(Fig.  16)  goes  to  zero  at  periods  of  low  tide  during  the  surge.  After  the  surge,  transport 
rates  for  breaches  at  z  =  1.5  m  and  z  =  2  m  continue  to  oscillate  at  the  tidal  frequency, 
indicating  an  open  breach  remains.  Intermittent  drying  of  the  inlet  causes  a  much 
reduced  transport  for  the  breach  channels  that  had  high  pilot  channels  (z  =  3  m,  z  = 
2.5  m)  where  little  to  no  inundation  is  present.  Because  of  the  reduced  transport,  the 
higher  elevation  channels  remain  closed  after  the  surge,  and  cross-sectional  area  is 
constant  (Fig.  17). 
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Fig.  14.  Calculated  sediment  transport  rate  through  the  four  breaches,  Case  Surge 


Fig.  15.  Breach  channel  cross-sectional  area  evolution  ,  Case  Surge 
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Fig.  16.  Calculated  sediment  transport  rate  through  the  four  breaches,  Case  Surge+M2 
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Fig.  17.  Breach  channel  cross-sectional  area  evolution  ,  Case  Surge+M2 


CONCLUDING  DISCUSSION 

This  study  investigated  predictions  of  a  new  numerical  model  of  multiple  coastal  barrier 
breaches  at  a  regional  scale  allowing,  for  example,  investigation  of  processes  associated 
with  independent  evolution  of  each  breach.  Numerical  simulations  of  coastal  breaches  at 
regional  scale  improve  understanding  of  the  physical  processes  by  isolating  forcing 
mechanisms  and  comparing  responses  of  breach  channels  of  varying  dimensions. 
Enhancements  to  the  morphologic  breaching  model  add  the  capability  to  simulate 
multiple  breach  development  at  the  regional  scale  with  boundary  condition  forcing 
interaction  with  Cascade. 
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Calculation  results  indicate  that  slow,  shallow  inundation  of  breach  pilot  channels  tends 
to  produce  breach  width  expanding  at  a  greater  rate  than  breach  depth.  This  is  likely  due 
to  the  increased  shear  stress  applied  over  a  wider  breach  channel  surface  area  in  the 
model.  A  shallow  breach  also  promotes  a  longer  incipient  breach  phase  and  longer  time 
before  developing  a  tidally  efficient  channel.  Fast  and  large  inundation  appears  to 
promote  rapid  increase  in  depth  as  well  as  a  more  rapid  incipient  breach  phase,  where  the 
cross-sectional  area  quickly  changes. 

Future  research  will  examine  model  performance  by  comparison  of  predictions  with 
observations  of  multiple  breaches  serving  the  same  bay.  Developments  of  the  regional 
breaching  model  are  expected  to  continue  with  an  embedded  form  in  Cascade  for  fully 
automated  forcing  and  interaction.  Finally,  interaction  between  breaches  serving  the 
same  bay  will  need  to  be  represented  for  examining  multiple  inlet  stability  for  a  nearly 
fixed  available  tidal  prism. 
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